# Final

library(readr)
## Warning: package 'readr' was built under R version 4.1.2
library(fpp)
## Loading required package: forecast
## Warning: package 'forecast' was built under R version 4.1.2
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
## Loading required package: fma
## Loading required package: expsmooth
## Loading required package: lmtest
## Warning: package 'lmtest' was built under R version 4.1.2
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.1.2
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## Loading required package: tseries
## Warning: package 'tseries' was built under R version 4.1.2
library(fpp2)
## ── Attaching packages ────────────────────────────────────────────── fpp2 2.4 ──
## ✔ ggplot2 3.3.5
## 
## 
## Attaching package: 'fpp2'
## The following objects are masked from 'package:fpp':
## 
##     ausair, ausbeer, austa, austourists, debitcards, departures,
##     elecequip, euretail, guinearice, oil, sunspotarea, usmelec
IPG3113N <- read_csv("/Users/sarthakmehta/Desktop/IPG3113N.csv")
## Rows: 121 Columns: 2
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## dbl  (1): IPG3113N
## date (1): DATE
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
candy_ts <- ts(IPG3113N$IPG3113N,frequency = 12,start=c(2012,10))

#Plot and Inference
candy_ts
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 2012                                                                        
## 2013  92.8393  88.4378  88.9572  82.6418  79.6009  79.7216  80.3623  87.3990
## 2014  91.8045  91.3464  88.8527  81.7442  77.6292  78.3359  77.0192  85.4452
## 2015  93.6676  92.4733  90.2388  85.5111  81.8907  85.4442  87.4143  99.0757
## 2016  97.1509  97.9475  98.4111  95.0657  93.3682  93.8015  94.3271  96.2487
## 2017 100.4828 103.5407  95.2095  96.9654  91.2276  92.9268  90.3917 100.6988
## 2018  99.9651 101.7932  97.6257  97.1776  93.2973  90.7030  90.7825  96.9555
## 2019 102.3790 103.2079 102.5133  96.6862  93.2061  93.6885  91.3219  99.4256
## 2020 100.1032 102.5297  96.9412  82.4777  84.4155  87.5411  91.2016 100.3417
## 2021 110.1841 101.7558 102.2030  94.1243  94.0966  94.4684  92.0391  98.7301
## 2022 113.1828 111.8384 110.4263 111.0788 103.3188 105.4497 104.2348 105.6166
##           Sep      Oct      Nov      Dec
## 2012           99.8541 100.7874 100.7643
## 2013  92.0523 102.5585 106.4783 108.0668
## 2014  92.7103 103.4097 110.9082 109.9962
## 2015  99.8525 110.3226 109.3626 106.4516
## 2016 101.3271 111.2323 108.1348 107.5859
## 2017 106.6639 106.8689 105.1628 109.8611
## 2018 105.8597 105.4468 106.0408 111.4674
## 2019 108.2320 109.2961 112.5268 115.7811
## 2020 106.3095 114.2083 116.2062 117.5581
## 2021 106.2289 126.5753 127.6255 129.0640
## 2022 119.0161 130.2894
plot(candy_ts) #plotting the candy_ts data

#Observation: Through the graph we can that the monthly production of candies in the US has been pretty cyclical means that we 
#can see a clear pattern to it. And there are sharp declines and rises throughout the data, so to cut out a window to forecast 
#would not be correct. Hence would be going on with the full data set provided.

#Central Tendency

min(candy_ts)
## [1] 77.0192
# [1] 77.0192

max(candy_ts)
## [1] 130.2894
# [1] 130.2894

mean(candy_ts)
## [1] 99.44454
# [1] 99.44454

median(candy_ts)
## [1] 99.8525
# [1] 99.8525

quantile(candy_ts)
##       0%      25%      50%      75%     100% 
##  77.0192  92.0523  99.8525 106.4516 130.2894
#   25%      75%    
# 92.0523  106.4516  

boxplot((candy_ts))

#Summary: Well the summaries show us that the amongst the median values calculated min value of candy production in the US
#period of 2012 to 2022 was 77.0192 and the max value was 130.2894. Furthermore, the mean value was 99.44454$ and 
#the median value was 99.8525. The 1st quartile was 92.0523 and the 3rd quartile was 106.4516. The box plot tells us that
#that the there are 2 outliers, and the minimum value is a bit more closer to the 1st quartile than the maximum value is to 
#the 3rd quatile,and the median is approximate in between the 1st and 3rd quartile. But when we look closely at the mean and 
#median we can see that as mean is a bit lower than the median (2nd quartile line) hence would me that data is negatively skewed.


#Decomposition

stl_decomp <- stl(candy_ts,s.window =12)
stl_decomp
##  Call:
##  stl(x = candy_ts, s.window = 12)
## 
## Components
##             seasonal     trend     remainder
## Oct 2012  10.5991242  88.76885   0.486128404
## Nov 2012  11.9082450  88.94684  -0.067689965
## Dec 2012  12.4396616  89.12484  -0.800204106
## Jan 2013   0.9586825  89.30284   2.577777391
## Feb 2013   0.5607168  89.50068  -1.623601801
## Mar 2013  -2.0268681  89.69853   1.285538278
## Apr 2013  -6.6117581  89.89637  -0.642816560
## May 2013 -10.0729539  90.10035  -0.426498736
## Jun 2013  -9.3367546  90.30433  -1.245976113
## Jul 2013  -9.5950782  90.50831  -0.550930554
## Aug 2013  -2.0496520  90.60483  -1.156182465
## Sep 2013   3.3932326  90.70136  -2.042292762
## Oct 2013  10.5381605  90.79789   1.222453628
## Nov 2013  11.7366414  90.68594   4.055718666
## Dec 2013  12.3879538  90.57399   5.104852222
## Jan 2014   0.9628192  90.46205   0.379632771
## Feb 2014   0.6389851  90.41520   0.292219599
## Mar 2014  -2.0111343  90.36834   0.495491725
## Apr 2014  -6.6100457  90.32149  -1.967244165
## May 2014 -10.0306151  90.46147  -2.801659713
## Jun 2014  -9.3114433  90.60146  -2.954116492
## Jul 2014  -9.5920828  90.74144  -4.130162024
## Aug 2014  -2.0542340  90.98619  -3.486757861
## Sep 2014   3.4924209  91.23094  -2.013059763
## Oct 2014  10.4813530  91.47569   1.452661089
## Nov 2014  11.5692114  91.99909   7.339900617
## Dec 2014  12.3404369  92.52249   5.133272951
## Jan 2015   0.9713533  93.04589  -0.349645527
## Feb 2015   0.7218572  93.63809  -1.886647844
## Mar 2015  -1.9905303  94.23029  -2.000958834
## Apr 2015  -6.6031967  94.82249  -2.708190937
## May 2015  -9.9830309  95.09530  -3.221568015
## Jun 2015  -9.2807779  95.36811  -0.643132352
## Jul 2015  -9.5840192  95.64092   1.357397749
## Aug 2015  -2.0540339  96.13102   4.998711111
## Sep 2015   3.5954289  96.62112  -0.364053081
## Oct 2015  10.4285374  97.11123   2.782837102
## Nov 2015  11.3025047  97.69738   0.362714208
## Dec 2015  12.2706801  98.28354  -4.102616766
## Jan 2016   0.9304921  98.86969  -2.649284402
## Feb 2016   0.8764821  99.08581  -2.014796319
## Mar 2016  -1.9993279  99.30194   1.108491712
## Apr 2016  -6.6292411  99.51806   2.176882941
## May 2016  -9.9089575  99.61689   3.660265343
## Jun 2016  -9.2001789  99.71573   3.285952692
## Jul 2016  -9.5406504  99.81456   4.053190219
## Aug 2016  -2.0359966  99.88701  -1.602309051
## Sep 2016   3.7863623  99.95945  -2.418713451
## Oct 2016  10.3551207 100.03190   0.845282712
## Nov 2016  11.0123956 100.00919  -2.886783341
## Dec 2016  12.1747194  99.98648  -4.575298361
## Jan 2017   0.8589920  99.96377  -0.339962196
## Feb 2017   0.9960331  99.98942   2.555249318
## Mar 2017  -2.0487293 100.01507  -2.756835789
## Apr 2017  -6.7014192 100.04071   3.626106613
## May 2017  -9.8860464 100.01204   1.101609349
## Jun 2017  -9.1757707  99.98336   2.119209238
## Jul 2017  -9.5560123  99.95469  -0.006973671
## Aug 2017  -2.0792298  99.95389   2.824143474
## Sep 2017   3.9180669  99.95309   2.792746411
## Oct 2017  10.5189014  99.95229  -3.602288426
## Nov 2017  11.0420904  99.94582  -5.825108760
## Dec 2017  12.4069792  99.93935  -2.485228857
## Jan 2018   1.0703083  99.93288  -1.038089254
## Feb 2018   0.9505060  99.87849   0.964202807
## Mar 2018  -2.0618904  99.82410  -0.136510979
## Apr 2018  -6.7759368  99.76971   4.183825153
## May 2018  -9.9575216  99.81734   3.437481672
## Jun 2018  -9.2362093  99.86497   0.074241035
## Jul 2018  -9.7580459  99.91260   0.627949363
## Aug 2018  -2.6602766 100.03749  -0.421708467
## Sep 2018   3.8733267 100.16237   1.823999816
## Oct 2018  10.7524738 100.28726  -5.592935852
## Nov 2018  11.1397392 100.43136  -5.530297287
## Dec 2018  12.7053552 100.57545  -1.813409278
## Jan 2019   1.3426650 100.71955   0.316784921
## Feb 2019   0.9609437 100.96702   1.279936357
## Mar 2019  -2.0251972 101.21449   3.324007384
## Apr 2019  -6.8067104 101.46196   2.030950738
## May 2019  -9.9908230 101.63893   1.557993363
## Jun 2019  -9.2640441 101.81590   1.136644525
## Jul 2019  -9.9317231 101.99287  -0.739246400
## Aug 2019  -3.2172141 101.63460   1.008214976
## Sep 2019   3.8494705 101.27633   3.106200634
## Oct 2019  10.9518486 100.91806  -2.573807143
## Nov 2019  11.2186645 100.45449   0.853641995
## Dec 2019  12.8851579  99.99093   2.905013508
## Jan 2020   1.4674880  99.52736  -0.891651539
## Feb 2020   0.9667147  99.49509   2.067894875
## Mar 2020  -1.9830794  99.46282  -0.538537889
## Apr 2020  -6.8148633  99.43054 -10.137980898
## May 2020 -10.0329043  99.74307  -5.294666481
## Jun 2020  -9.3135889 100.05560  -3.200908472
## Jul 2020 -10.0346542 100.36812   0.868130191
## Aug 2020  -3.4609873 100.86804   2.934642849
## Sep 2020   3.8269911 101.36796   1.114543976
## Oct 2020  11.1280938 101.86789   1.212320792
## Nov 2020  11.2760261 102.34620   2.583971116
## Dec 2020  13.0449631 102.82452   1.688616745
## Jan 2021   1.5755494 103.30284   5.305713092
## Feb 2021   0.9589601 103.67269  -2.875850112
## Mar 2021  -1.9503607 104.04254   0.110818098
## Apr 2021  -6.8282886 104.41240  -3.459806511
## May 2021 -10.0761009 105.08146  -0.908754977
## Jun 2021  -9.3600919 105.75052  -1.922024797
## Jul 2021 -10.1312918 106.41958  -4.249185771
## Aug 2021  -3.6952153 107.27536  -4.850048264
## Sep 2021   3.8158456 108.13115  -5.718095113
## Oct 2021  11.2130846 108.98694   6.375279931
## Nov 2021  11.3363315 109.81832   6.470852045
## Dec 2021  13.1803233 110.64970   5.233979330
## Jan 2022   1.6748211 111.48108   0.026900529
## Feb 2022   0.9604557 112.09593  -1.217984182
## Mar 2022  -1.9154152 112.71078  -0.369063354
## Apr 2022  -6.8967372 113.32563   4.649908547
## May 2022 -10.1274938 113.89493  -0.448631708
## Jun 2022  -9.4097104 114.46422   0.395188008
## Jul 2022 -10.2055011 115.03352  -0.593218213
## Aug 2022  -3.8402816 115.52604  -6.069163191
## Sep 2022   3.7907091 116.01857  -0.793179360
## Oct 2022  11.2961746 116.51110   2.482129583
plot(stl_decomp) #Decomposition Plot

#Yes the time series is seasonal, and it is additive seasonal not multiplicative as the magnitude is constant. Initially from 
#the month of October - December had pretty high seasonal values and later we also see that this trend is continuing and the 
#months of October - December have really high seasonal values. I believe the reason for it is that around that time of the year 
#we have the holidays season meaning thanksgiving, Christmas and new years so naturally the consumers would be needing more 
#candies to celebrate the holidays. As for the rest of the years as expected the seasonal values are negative because no season 
#is like the holiday season and the production during those months is really high therefore in the months following the 
#seasonal values tend to become very low or even negative.

seasadj(stl_decomp)
##            Jan       Feb       Mar       Apr       May       Jun       Jul
## 2012                                                                      
## 2013  91.88062  87.87708  90.98407  89.25356  89.67385  89.05835  89.95738
## 2014  90.84168  90.70741  90.86383  88.35425  87.65982  87.64734  86.61128
## 2015  92.69625  91.75144  92.22933  92.11430  91.87373  94.72498  96.99832
## 2016  96.22041  97.07102 100.41043 101.69494 103.27716 103.00168 103.86775
## 2017  99.62381 102.54467  97.25823 103.66682 101.11365 102.10257  99.94771
## 2018  98.89479 100.84269  99.68759 103.95354 103.25482  99.93921 100.54055
## 2019 101.03634 102.24696 104.53850 103.49291 103.19692 102.95254 101.25362
## 2020  98.63571 101.56299  98.92428  89.29256  94.44840  96.85469 101.23625
## 2021 108.60855 100.79684 104.15336 100.95259 104.17270 103.82849 102.17039
## 2022 111.50798 110.87794 112.34172 117.97554 113.44629 114.85941 114.44030
##            Aug       Sep       Oct       Nov       Dec
## 2012                      89.25498  88.87915  88.32464
## 2013  89.44865  88.65907  92.02034  94.74166  95.67885
## 2014  87.49943  89.21788  92.92835  99.33899  97.65576
## 2015 101.12973  96.25707  99.89406  98.06010  94.18092
## 2016  98.28470  97.54074 100.87718  97.12240  95.41118
## 2017 102.77803 102.74583  96.35000  94.12071  97.45412
## 2018  99.61578 101.98637  94.69433  94.90106  98.76204
## 2019 102.64281 104.38253  98.34425 101.30814 102.89594
## 2020 103.80269 102.48251 103.08021 104.93017 104.51314
## 2021 102.42532 102.41305 115.36222 116.28917 115.88368
## 2022 109.45688 115.22539 118.99323
plot(candy_ts)
lines(seasadj(stl_decomp), col="Orange")

#After taking out the seasonal component the graph still looks the same in the overall direction and trend but the values have
#shrunk dramatically this means that the magnitude of values in the timeseries is affected but the seasonal component. Yes, 
#seasonality has big fluctuations to the values of the timeseries.


#Naive

naive_forecast <- naive(candy_ts,12) #give me the forecast for the next 12 months using Naive forecast method
plot(naive_forecast)

attributes(naive_forecast)
## $names
##  [1] "method"    "model"     "lambda"    "x"         "fitted"    "residuals"
##  [7] "series"    "mean"      "level"     "lower"     "upper"    
## 
## $class
## [1] "forecast"
#$names
#[1] "method"    "model"     "lambda"    "x"         "fitted"    "residuals" "series"    "mean"      "level"     "lower"    
#[11] "upper"    

#$class
#[1] "forecast"

plot(naive_forecast$residuals) #the residual plot shows 6 anomalies(not within range from what I can see); 1st one being at the

#start of year 2014 where there is a sharp decline in production same goes for the start of year 2015,2020,2022 and in the 
#beginning months of year 2022 and the 6th one being the sharp increase in the production just at the end of year 2021.

hist(naive_forecast$residuals) #shows a pretty decent normal distribution

Acf(naive_forecast$residuals) #shows most values are out of the threshold region meaning they are dependent on the previous 

#values (correaltion present).

accuracy(naive_forecast)
##                     ME     RMSE      MAE        MPE     MAPE    MASE      ACF1
## Training set 0.2536275 6.305638 4.625101 0.02809525 4.635846 1.05288 0.2207691
#                 ME     RMSE      MAE        MPE      MAPE       MASE      ACF1
#Training set 0.2536275 6.305638 4.625101 0.02809525 4.635846   1.05288   0.2207691


ets_forecast <- ets(candy_ts)
plot(candy_ts)

forecast_ets <- forecast(ets_forecast, h=12) #forecasting for the next 12 months
plot(forecast_ets)

forecast_ets
##          Point Forecast     Lo 80    Hi 80     Lo 95    Hi 95
## Nov 2022       129.6266 123.89521 135.3579 120.86121 138.3919
## Dec 2022       131.4585 124.44955 138.4674 120.73925 142.1777
## Jan 2023       120.2379 112.53103 127.9448 108.45126 132.0245
## Feb 2023       119.7004 111.19327 128.2074 106.68990 132.7108
## Mar 2023       117.3170 108.13076 126.5032 103.26787 131.3661
## Apr 2023       112.4141 102.67755 122.1506  97.52333 127.3049
## May 2023       109.4166  99.16483 119.6684  93.73787 125.0953
## Jun 2023       110.4319  99.63670 121.2271  93.92206 126.9418
## Jul 2023       110.2762  98.97651 121.5760  92.99479 127.5577
## Aug 2023       117.4103 105.49177 129.3289  99.18246 135.6382
## Sep 2023       123.6545 111.09743 136.2116 104.45011 142.8589
## Oct 2023       131.8069 118.54708 145.0667 111.52776 152.0860
#        Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
#Nov 2022       129.6266 123.89521 135.3579 120.86121 138.3919
#Dec 2022       131.4585 124.44955 138.4674 120.73925 142.1777
#Jan 2023       120.2379 112.53103 127.9448 108.45126 132.0245
#Feb 2023       119.7004 111.19327 128.2074 106.68990 132.7108
#Mar 2023       117.3170 108.13076 126.5032 103.26787 131.3661
#Apr 2023       112.4141 102.67755 122.1506  97.52333 127.3049
#May 2023       109.4166  99.16483 119.6684  93.73787 125.0953
#Jun 2023       110.4319  99.63670 121.2271  93.92206 126.9418
#Jul 2023       110.2762  98.97651 121.5760  92.99479 127.5577
#Aug 2023       117.4103 105.49177 129.3289  99.18246 135.6382
#Sep 2023       123.6545 111.09743 136.2116 104.45011 142.8589
#Oct 2023       131.8069 118.54708 145.0667 111.52776 152.0860

ets_forecast$mse
## [1] 10.95792
#[1] 10.95792
#RSME is 6.305638 on average forecast values were 6.305638 away from the actual
#MPE is 2.8% percenatge of error is around 2.8%


#Simple Moving Average

plot(candy_ts)
MA3_forecast <- ma(candy_ts,order=3) #taking into account the 5 most recent values
MA6_forecast <- ma(candy_ts,order=6) #taking into account the 9 most recent values
MA9_forecast <- ma(candy_ts,order=9) #taking into account the 9 most recent values


rwf_forecast <- rwf(candy_ts,12)
snaive_forecast <- snaive(candy_ts, 12)

plot(candy_ts)
lines(rwf_forecast$mean,col="Orange")
lines(snaive_forecast$mean,col="brown") #works best for my time series: snaive method
lines(MA3_forecast,col="Red")
lines(MA6_forecast,col="Blue")
lines(MA9_forecast,col="Green")

#as the order goes up the line comes up shorter and shorter from the end, and from the start as well, the sample in getting 
#smaller. And also as the order goes up the trend is somewhat followed but the magnitude shrinks in size as well.

#Simple Smoothing

SSE_Simple <- holt(candy_ts,gamma=FALSE)
attributes(SSE_Simple)
## $names
##  [1] "model"     "mean"      "level"     "x"         "upper"     "lower"    
##  [7] "fitted"    "method"    "series"    "residuals"
## 
## $class
## [1] "forecast"
plot(SSE_Simple)

attributes(SSE_Simple)
## $names
##  [1] "model"     "mean"      "level"     "x"         "upper"     "lower"    
##  [7] "fitted"    "method"    "series"    "residuals"
## 
## $class
## [1] "forecast"
#$names
#[1] "model"     "mean"      "level"     "x"         "upper"     "lower"     "fitted"    "method"    "series"    "residuals"

#$class
#[1] "forecast"

Acf(SSE_Simple$residuals)

Box.test(residuals(SSE_Simple), lag=12, type="Ljung")
## 
##  Box-Ljung test
## 
## data:  residuals(SSE_Simple)
## X-squared = 147.8, df = 12, p-value < 2.2e-16
#data:  residuals(auto_fit)
#X-squared = 13.084, df = 12, p-value < 2.2e-16

#As we can see from the Ljung box statistic that the p-values are < 2.2e-16 which makes them significant and we accept the null 
#hypothesis which means that the residual values are dependent

plot.ts(residuals(SSE_Simple))

hist(SSE_Simple$residuals)  #normal distribution

accuracy(SSE_Simple)
##                      ME    RMSE      MAE        MPE     MAPE     MASE      ACF1
## Training set 0.04457017 6.29552 4.615979 -0.1819999 4.631386 1.050803 0.2185404
#                 ME     RMSE      MAE        MPE     MAPE     MASE      ACF1
#Training set 0.04457017 6.29552 4.615979 -0.1819999 4.631386 1.050803 0.2185404

SSE_Simple$mse
## NULL
#[1] NULL
#RSME is 6.29552 on average forecast values were 6.29552 away from the actual
#MPE is 18% percenatge of error is around 18%


#Holt Winters

HW_forecast <- HoltWinters(candy_ts)
plot(HW_forecast)

HW_forecast
## Holt-Winters exponential smoothing with trend and additive seasonal component.
## 
## Call:
## HoltWinters(x = candy_ts)
## 
## Smoothing parameters:
##  alpha: 0.7137653
##  beta : 0
##  gamma: 0.3677262
## 
## Coefficients:
##             [,1]
## a   118.58175151
## b     0.02735169
## s1   12.73645821
## s2   15.64347958
## s3    2.86140879
## s4    1.55347367
## s5   -1.22825112
## s6   -6.63034922
## s7  -10.18821076
## s8   -9.91381693
## s9  -10.71852359
## s10  -3.99912312
## s11   2.78200916
## s12  10.83734441
HWForecast <- predict(HW_forecast, 12)
plot(HWForecast)

HWForecast
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 2022                                                                        
## 2023 121.5252 120.2446 117.4903 112.1155 108.5850 108.8867 108.1094 114.8561
##           Sep      Oct      Nov      Dec
## 2022                   131.3456 134.2799
## 2023 121.6646 129.7473
HW_forecast #Alpha = 0.7137653 (the newest values does not have maximum weight but has relatively more weight than the last value)
## Holt-Winters exponential smoothing with trend and additive seasonal component.
## 
## Call:
## HoltWinters(x = candy_ts)
## 
## Smoothing parameters:
##  alpha: 0.7137653
##  beta : 0
##  gamma: 0.3677262
## 
## Coefficients:
##             [,1]
## a   118.58175151
## b     0.02735169
## s1   12.73645821
## s2   15.64347958
## s3    2.86140879
## s4    1.55347367
## s5   -1.22825112
## s6   -6.63034922
## s7  -10.18821076
## s8   -9.91381693
## s9  -10.71852359
## s10  -3.99912312
## s11   2.78200916
## s12  10.83734441
#beta : 0 (this is the coefficient of trend smoothing in HW)
#gamma: 0.3677262 (seasonal model)

ets_forecast #sigma: 0.0345 (standard dev. of residuals)
## ETS(M,A,A) 
## 
## Call:
##  ets(y = candy_ts) 
## 
##   Smoothing parameters:
##     alpha = 0.6829 
##     beta  = 1e-04 
##     gamma = 1e-04 
## 
##   Initial states:
##     l = 90.6178 
##     b = 0.2533 
##     s = 3.0373 -2.9535 -9.8347 -9.4264 -10.1889 -6.9389
##            -1.7831 0.8532 1.6436 13.117 11.537 10.9375
## 
##   sigma:  0.0345
## 
##      AIC     AICc      BIC 
## 894.2047 900.1464 941.7331
#Initial states:l = 90.6178  
#b = 0.2533 

# AIC     AICc      BIC 
#894.2047 900.1464 941.7331 

attributes(HW_forecast)
## $names
## [1] "fitted"       "x"            "alpha"        "beta"         "gamma"       
## [6] "coefficients" "seasonal"     "SSE"          "call"        
## 
## $class
## [1] "HoltWinters"
#$names
#[1] "fitted"       "x"            "alpha"        "beta"         "gamma"        "coefficients" "seasonal"     "SSE"         
#[9] "call"        

#$class
#[1] "HoltWinters"

attributes(HWForecast)
## $dim
## [1] 12  1
## 
## $dimnames
## $dimnames[[1]]
## NULL
## 
## $dimnames[[2]]
## [1] "fit"
## 
## 
## $tsp
## [1] 2022.833 2023.750   12.000
## 
## $class
## [1] "ts"
#$dim
#[1] 12  1

#$dimnames
#$dimnames[[1]]
#NULL

#$dimnames[[2]]
#[1] "fit"

#$tsp
#[1] 2022.833 2023.750   12.000

#$class
#[1] "ts"

HWFResiduals <-residuals(HW_forecast)

Acf(HWFResiduals) #most values seem to be in the threshold region meaning apart from the 2 values outside of the region all 

#values are dependent on the last value

Box.test(HWFResiduals, lag=12, type="Ljung")
## 
##  Box-Ljung test
## 
## data:  HWFResiduals
## X-squared = 27.423, df = 12, p-value = 0.006713
#Box-Ljung test
#data:  HWFResiduals
#X-squared = 27.423, df = 12, p-value = 0.006713

#As we can see from the Ljung box statistic that the p-values are 0.006713 which makes them significant and we accept the null 
#hypothesis which means that the residual values are dependent.

plot.ts(HWFResiduals)

hist(HWFResiduals)  #normal distribution

#ARIMA or Box-Jenkins

plot(candy_ts)

ndiffs(candy_ts) #as the ndiffs function gives out the value 1 it means that the model is not stationary, and that would mean 
## [1] 1
#that we need 1 difference to get the stationary value.
#Seasonality is not needed for this.

tsdisplay(candy_ts)

candydiff1 <- diff(candy_ts, differences=1)
plot(candydiff1)

tsdisplay(candydiff1)

auto_fit <- auto.arima(candy_ts, trace=TRUE, stepwise = FALSE)
## 
##  ARIMA(0,0,0)(0,1,0)[12]                    : 695.7015
##  ARIMA(0,0,0)(0,1,0)[12] with drift         : 670.7642
##  ARIMA(0,0,0)(0,1,1)[12]                    : 695.999
##  ARIMA(0,0,0)(0,1,1)[12] with drift         : Inf
##  ARIMA(0,0,0)(0,1,2)[12]                    : 696.3967
##  ARIMA(0,0,0)(0,1,2)[12] with drift         : Inf
##  ARIMA(0,0,0)(1,1,0)[12]                    : 696.3465
##  ARIMA(0,0,0)(1,1,0)[12] with drift         : 672.2067
##  ARIMA(0,0,0)(1,1,1)[12]                    : Inf
##  ARIMA(0,0,0)(1,1,1)[12] with drift         : Inf
##  ARIMA(0,0,0)(1,1,2)[12]                    : Inf
##  ARIMA(0,0,0)(1,1,2)[12] with drift         : Inf
##  ARIMA(0,0,0)(2,1,0)[12]                    : 697.5285
##  ARIMA(0,0,0)(2,1,0)[12] with drift         : 666.0752
##  ARIMA(0,0,0)(2,1,1)[12]                    : Inf
##  ARIMA(0,0,0)(2,1,1)[12] with drift         : Inf
##  ARIMA(0,0,0)(2,1,2)[12]                    : Inf
##  ARIMA(0,0,0)(2,1,2)[12] with drift         : Inf
##  ARIMA(0,0,1)(0,1,0)[12]                    : 654.6188
##  ARIMA(0,0,1)(0,1,0)[12] with drift         : 639.2968
##  ARIMA(0,0,1)(0,1,1)[12]                    : 656.1893
##  ARIMA(0,0,1)(0,1,1)[12] with drift         : Inf
##  ARIMA(0,0,1)(0,1,2)[12]                    : 657.6255
##  ARIMA(0,0,1)(0,1,2)[12] with drift         : Inf
##  ARIMA(0,0,1)(1,1,0)[12]                    : 656.284
##  ARIMA(0,0,1)(1,1,0)[12] with drift         : 637.5419
##  ARIMA(0,0,1)(1,1,1)[12]                    : 657.9208
##  ARIMA(0,0,1)(1,1,1)[12] with drift         : Inf
##  ARIMA(0,0,1)(1,1,2)[12]                    : 659.7173
##  ARIMA(0,0,1)(1,1,2)[12] with drift         : Inf
##  ARIMA(0,0,1)(2,1,0)[12]                    : 657.8574
##  ARIMA(0,0,1)(2,1,0)[12] with drift         : 634.4358
##  ARIMA(0,0,1)(2,1,1)[12]                    : 659.7902
##  ARIMA(0,0,1)(2,1,1)[12] with drift         : Inf
##  ARIMA(0,0,1)(2,1,2)[12]                    : Inf
##  ARIMA(0,0,1)(2,1,2)[12] with drift         : Inf
##  ARIMA(0,0,2)(0,1,0)[12]                    : 636.3448
##  ARIMA(0,0,2)(0,1,0)[12] with drift         : 626.0602
##  ARIMA(0,0,2)(0,1,1)[12]                    : 637.898
##  ARIMA(0,0,2)(0,1,1)[12] with drift         : Inf
##  ARIMA(0,0,2)(0,1,2)[12]                    : 636.7653
##  ARIMA(0,0,2)(0,1,2)[12] with drift         : Inf
##  ARIMA(0,0,2)(1,1,0)[12]                    : 638.0881
##  ARIMA(0,0,2)(1,1,0)[12] with drift         : 625.6669
##  ARIMA(0,0,2)(1,1,1)[12]                    : 638.334
##  ARIMA(0,0,2)(1,1,1)[12] with drift         : Inf
##  ARIMA(0,0,2)(1,1,2)[12]                    : 638.992
##  ARIMA(0,0,2)(1,1,2)[12] with drift         : Inf
##  ARIMA(0,0,2)(2,1,0)[12]                    : 637.8508
##  ARIMA(0,0,2)(2,1,0)[12] with drift         : 621.3809
##  ARIMA(0,0,2)(2,1,1)[12]                    : 638.4469
##  ARIMA(0,0,2)(2,1,1)[12] with drift         : Inf
##  ARIMA(0,0,3)(0,1,0)[12]                    : 629.9275
##  ARIMA(0,0,3)(0,1,0)[12] with drift         : 622.7562
##  ARIMA(0,0,3)(0,1,1)[12]                    : 629.7434
##  ARIMA(0,0,3)(0,1,1)[12] with drift         : Inf
##  ARIMA(0,0,3)(0,1,2)[12]                    : 629.0504
##  ARIMA(0,0,3)(0,1,2)[12] with drift         : Inf
##  ARIMA(0,0,3)(1,1,0)[12]                    : 630.4963
##  ARIMA(0,0,3)(1,1,0)[12] with drift         : 620.9963
##  ARIMA(0,0,3)(1,1,1)[12]                    : 629.6219
##  ARIMA(0,0,3)(1,1,1)[12] with drift         : Inf
##  ARIMA(0,0,3)(2,1,0)[12]                    : 630.7565
##  ARIMA(0,0,3)(2,1,0)[12] with drift         : 618.3351
##  ARIMA(0,0,4)(0,1,0)[12]                    : 629.5131
##  ARIMA(0,0,4)(0,1,0)[12] with drift         : 623.8854
##  ARIMA(0,0,4)(0,1,1)[12]                    : 628.2313
##  ARIMA(0,0,4)(0,1,1)[12] with drift         : Inf
##  ARIMA(0,0,4)(1,1,0)[12]                    : 629.5068
##  ARIMA(0,0,4)(1,1,0)[12] with drift         : 621.9546
##  ARIMA(0,0,5)(0,1,0)[12]                    : 630.6005
##  ARIMA(0,0,5)(0,1,0)[12] with drift         : 625.8286
##  ARIMA(1,0,0)(0,1,0)[12]                    : 622.9559
##  ARIMA(1,0,0)(0,1,0)[12] with drift         : 619.2546
##  ARIMA(1,0,0)(0,1,1)[12]                    : 616.096
##  ARIMA(1,0,0)(0,1,1)[12] with drift         : Inf
##  ARIMA(1,0,0)(0,1,2)[12]                    : Inf
##  ARIMA(1,0,0)(0,1,2)[12] with drift         : Inf
##  ARIMA(1,0,0)(1,1,0)[12]                    : 620.045
##  ARIMA(1,0,0)(1,1,0)[12] with drift         : 615.4554
##  ARIMA(1,0,0)(1,1,1)[12]                    : Inf
##  ARIMA(1,0,0)(1,1,1)[12] with drift         : Inf
##  ARIMA(1,0,0)(1,1,2)[12]                    : Inf
##  ARIMA(1,0,0)(1,1,2)[12] with drift         : Inf
##  ARIMA(1,0,0)(2,1,0)[12]                    : 619.4956
##  ARIMA(1,0,0)(2,1,0)[12] with drift         : 613.3055
##  ARIMA(1,0,0)(2,1,1)[12]                    : Inf
##  ARIMA(1,0,0)(2,1,1)[12] with drift         : Inf
##  ARIMA(1,0,0)(2,1,2)[12]                    : Inf
##  ARIMA(1,0,0)(2,1,2)[12] with drift         : Inf
##  ARIMA(1,0,1)(0,1,0)[12]                    : 621.839
##  ARIMA(1,0,1)(0,1,0)[12] with drift         : 619.5991
##  ARIMA(1,0,1)(0,1,1)[12]                    : 615.267
##  ARIMA(1,0,1)(0,1,1)[12] with drift         : Inf
##  ARIMA(1,0,1)(0,1,2)[12]                    : Inf
##  ARIMA(1,0,1)(0,1,2)[12] with drift         : Inf
##  ARIMA(1,0,1)(1,1,0)[12]                    : 618.9704
##  ARIMA(1,0,1)(1,1,0)[12] with drift         : 616.1254
##  ARIMA(1,0,1)(1,1,1)[12]                    : Inf
##  ARIMA(1,0,1)(1,1,1)[12] with drift         : Inf
##  ARIMA(1,0,1)(1,1,2)[12]                    : Inf
##  ARIMA(1,0,1)(1,1,2)[12] with drift         : Inf
##  ARIMA(1,0,1)(2,1,0)[12]                    : 618.2582
##  ARIMA(1,0,1)(2,1,0)[12] with drift         : 614.2405
##  ARIMA(1,0,1)(2,1,1)[12]                    : Inf
##  ARIMA(1,0,1)(2,1,1)[12] with drift         : Inf
##  ARIMA(1,0,2)(0,1,0)[12]                    : 623.7407
##  ARIMA(1,0,2)(0,1,0)[12] with drift         : 620.8486
##  ARIMA(1,0,2)(0,1,1)[12]                    : Inf
##  ARIMA(1,0,2)(0,1,1)[12] with drift         : Inf
##  ARIMA(1,0,2)(0,1,2)[12]                    : Inf
##  ARIMA(1,0,2)(0,1,2)[12] with drift         : Inf
##  ARIMA(1,0,2)(1,1,0)[12]                    : 621.1253
##  ARIMA(1,0,2)(1,1,0)[12] with drift         : 618.2083
##  ARIMA(1,0,2)(1,1,1)[12]                    : Inf
##  ARIMA(1,0,2)(1,1,1)[12] with drift         : Inf
##  ARIMA(1,0,2)(2,1,0)[12]                    : 620.3847
##  ARIMA(1,0,2)(2,1,0)[12] with drift         : 616.1228
##  ARIMA(1,0,3)(0,1,0)[12]                    : 624.0376
##  ARIMA(1,0,3)(0,1,0)[12] with drift         : 622.9934
##  ARIMA(1,0,3)(0,1,1)[12]                    : Inf
##  ARIMA(1,0,3)(0,1,1)[12] with drift         : Inf
##  ARIMA(1,0,3)(1,1,0)[12]                    : 621.5939
##  ARIMA(1,0,3)(1,1,0)[12] with drift         : 620.4402
##  ARIMA(1,0,4)(0,1,0)[12]                    : 625.2055
##  ARIMA(1,0,4)(0,1,0)[12] with drift         : 624.5141
##  ARIMA(2,0,0)(0,1,0)[12]                    : 621.502
##  ARIMA(2,0,0)(0,1,0)[12] with drift         : 619.1749
##  ARIMA(2,0,0)(0,1,1)[12]                    : 615.6749
##  ARIMA(2,0,0)(0,1,1)[12] with drift         : Inf
##  ARIMA(2,0,0)(0,1,2)[12]                    : Inf
##  ARIMA(2,0,0)(0,1,2)[12] with drift         : Inf
##  ARIMA(2,0,0)(1,1,0)[12]                    : 619.0898
##  ARIMA(2,0,0)(1,1,0)[12] with drift         : 616.0103
##  ARIMA(2,0,0)(1,1,1)[12]                    : Inf
##  ARIMA(2,0,0)(1,1,1)[12] with drift         : Inf
##  ARIMA(2,0,0)(1,1,2)[12]                    : Inf
##  ARIMA(2,0,0)(1,1,2)[12] with drift         : Inf
##  ARIMA(2,0,0)(2,1,0)[12]                    : 618.3519
##  ARIMA(2,0,0)(2,1,0)[12] with drift         : 614.0563
##  ARIMA(2,0,0)(2,1,1)[12]                    : Inf
##  ARIMA(2,0,0)(2,1,1)[12] with drift         : Inf
##  ARIMA(2,0,1)(0,1,0)[12]                    : 623.5522
##  ARIMA(2,0,1)(0,1,0)[12] with drift         : 621.0785
##  ARIMA(2,0,1)(0,1,1)[12]                    : Inf
##  ARIMA(2,0,1)(0,1,1)[12] with drift         : Inf
##  ARIMA(2,0,1)(0,1,2)[12]                    : Inf
##  ARIMA(2,0,1)(0,1,2)[12] with drift         : Inf
##  ARIMA(2,0,1)(1,1,0)[12]                    : 619.2936
##  ARIMA(2,0,1)(1,1,0)[12] with drift         : 618.2342
##  ARIMA(2,0,1)(1,1,1)[12]                    : Inf
##  ARIMA(2,0,1)(1,1,1)[12] with drift         : Inf
##  ARIMA(2,0,1)(2,1,0)[12]                    : 617.3422
##  ARIMA(2,0,1)(2,1,0)[12] with drift         : Inf
##  ARIMA(2,0,2)(0,1,0)[12]                    : 625.7487
##  ARIMA(2,0,2)(0,1,0)[12] with drift         : 623.077
##  ARIMA(2,0,2)(0,1,1)[12]                    : Inf
##  ARIMA(2,0,2)(0,1,1)[12] with drift         : Inf
##  ARIMA(2,0,2)(1,1,0)[12]                    : 621.2392
##  ARIMA(2,0,2)(1,1,0)[12] with drift         : 620.4897
##  ARIMA(2,0,3)(0,1,0)[12]                    : 625.4338
##  ARIMA(2,0,3)(0,1,0)[12] with drift         : 624.7418
##  ARIMA(3,0,0)(0,1,0)[12]                    : 623.565
##  ARIMA(3,0,0)(0,1,0)[12] with drift         : 620.9597
##  ARIMA(3,0,0)(0,1,1)[12]                    : 617.5374
##  ARIMA(3,0,0)(0,1,1)[12] with drift         : Inf
##  ARIMA(3,0,0)(0,1,2)[12]                    : Inf
##  ARIMA(3,0,0)(0,1,2)[12] with drift         : Inf
##  ARIMA(3,0,0)(1,1,0)[12]                    : 621.2407
##  ARIMA(3,0,0)(1,1,0)[12] with drift         : 618.2282
##  ARIMA(3,0,0)(1,1,1)[12]                    : Inf
##  ARIMA(3,0,0)(1,1,1)[12] with drift         : Inf
##  ARIMA(3,0,0)(2,1,0)[12]                    : Inf
##  ARIMA(3,0,0)(2,1,0)[12] with drift         : 616.1571
##  ARIMA(3,0,1)(0,1,0)[12]                    : Inf
##  ARIMA(3,0,1)(0,1,0)[12] with drift         : Inf
##  ARIMA(3,0,1)(0,1,1)[12]                    : Inf
##  ARIMA(3,0,1)(0,1,1)[12] with drift         : Inf
##  ARIMA(3,0,1)(1,1,0)[12]                    : Inf
##  ARIMA(3,0,1)(1,1,0)[12] with drift         : 620.2729
##  ARIMA(3,0,2)(0,1,0)[12]                    : Inf
##  ARIMA(3,0,2)(0,1,0)[12] with drift         : Inf
##  ARIMA(4,0,0)(0,1,0)[12]                    : 625.6611
##  ARIMA(4,0,0)(0,1,0)[12] with drift         : 623.1987
##  ARIMA(4,0,0)(0,1,1)[12]                    : 619.4666
##  ARIMA(4,0,0)(0,1,1)[12] with drift         : Inf
##  ARIMA(4,0,0)(1,1,0)[12]                    : 623.4068
##  ARIMA(4,0,0)(1,1,0)[12] with drift         : 620.4931
##  ARIMA(4,0,1)(0,1,0)[12]                    : Inf
##  ARIMA(4,0,1)(0,1,0)[12] with drift         : Inf
##  ARIMA(5,0,0)(0,1,0)[12]                    : 626.6308
##  ARIMA(5,0,0)(0,1,0)[12] with drift         : 624.8284
## 
## 
## 
##  Best model: ARIMA(1,0,0)(2,1,0)[12] with drift
# Best model: ARIMA(1,0,0)(2,1,0)[12] with drift

auto_fit
## Series: candy_ts 
## ARIMA(1,0,0)(2,1,0)[12] with drift 
## 
## Coefficients:
##          ar1     sar1     sar2   drift
##       0.6441  -0.3042  -0.2351  0.2083
## s.e.  0.0742   0.1040   0.1099  0.0579
## 
## sigma^2 = 14.95:  log likelihood = -301.36
## AIC=612.72   AICc=613.31   BIC=626.18
# Series: candy_ts 
#ARIMA(1,0,0)(2,1,0)[12] with drift 

#Coefficients:
#        ar1     sar1     sar2   drift
#      0.6441  -0.3042  -0.2351  0.2083
#s.e.  0.0742   0.1040   0.1099  0.0579

#sigma^2 = 14.95:  log likelihood = -301.36
#AIC=612.72   AICc=613.31   BIC=626.18


attributes(auto_fit)
## $names
##  [1] "coef"      "sigma2"    "var.coef"  "mask"      "loglik"    "aic"      
##  [7] "arma"      "residuals" "call"      "series"    "code"      "n.cond"   
## [13] "nobs"      "model"     "xreg"      "bic"       "aicc"      "x"        
## [19] "fitted"   
## 
## $class
## [1] "forecast_ARIMA" "ARIMA"          "Arima"
#$names
#[1] "coef"      "sigma2"    "var.coef"  "mask"      "loglik"    "aic"       "arma"      "residuals" "call"      "series"   
#[11] "code"      "n.cond"    "nobs"      "model"     "xreg"      "bic"       "aicc"      "x"         "fitted"   

#$class
#[1] "forecast_ARIMA" "ARIMA"          "Arima" 


Acf(auto_fit$residuals)

Box.test(residuals(auto_fit), lag=12, type="Ljung")
## 
##  Box-Ljung test
## 
## data:  residuals(auto_fit)
## X-squared = 13.084, df = 12, p-value = 0.363
#data:  residuals(auto_fit)
#X-squared = 13.084, df = 12, p-value = 0.363

#As we can see from the Ljung box statistic that the p-values are 0.363 which makes them non-significant and we reject the null 
#hypothesis which means that the residual values are independent

plot.ts(residuals(auto_fit))

hist(auto_fit$residuals)  #normal distribution

tsdiag(auto_fit) #ACF test shows no relative significance amongst residuals as well same as the p-values for the Ljung-Box

accuracy(auto_fit)
##                        ME     RMSE      MAE        MPE     MAPE      MASE
## Training set -0.004444951 3.602064 2.583457 -0.1100037 2.571757 0.5881103
##                     ACF1
## Training set -0.06809078
#                   ME        RMSE      MAE      MPE      MAPE      MASE        ACF1
#Training set -0.004444951 3.602064 2.583457 -0.1100037 2.571757 0.5881103 -0.06809078

ARIMAF12 <- plot(forecast(auto_fit,h=12,level=c(99.5)))

ARIMAF12
## $mean
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 2022                                                                        
## 2023 115.0259 113.6240 111.0655 107.3725 102.3049 104.4700 104.2670 107.8069
##           Sep      Oct      Nov      Dec
## 2022                   130.2152 130.9781
## 2023 119.0309 130.1239                  
## 
## $lower
##            Jan       Feb       Mar       Apr       May       Jun       Jul
## 2022                                                                      
## 2023 101.35264  99.64659  96.96380  93.21964  88.13085  90.28708  90.08048
##            Aug       Sep       Oct       Nov       Dec
## 2022                               119.36102 118.06743
## 2023  93.61892 104.84220 115.93496                    
## 
## $upper
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 2022                                                                        
## 2023 128.6992 127.6015 125.1672 121.5254 116.4790 118.6528 118.4535 121.9950
##           Sep      Oct      Nov      Dec
## 2022                   141.0694 143.8888
## 2023 133.2195 144.3128
ARIMAF24 <-plot(forecast(auto_fit,h=24,level=c(99.5)))

ARIMAF24
## $mean
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 2022                                                                        
## 2023 115.0259 113.6240 111.0655 107.3725 102.3049 104.4700 104.2670 107.8069
## 2024 117.6141 114.5620 112.7877 108.3628 104.2936 106.0342 105.2377 109.3693
##           Sep      Oct      Nov      Dec
## 2022                   130.2152 130.9781
## 2023 119.0309 130.1239 130.6058 131.5482
## 2024 119.8675 133.1484                  
## 
## $lower
##            Jan       Feb       Mar       Apr       May       Jun       Jul
## 2022                                                                      
## 2023 101.35264  99.64659  96.96380  93.21964  88.13085  90.28708  90.08048
## 2024 100.49160  97.31920  95.49535  91.04988  86.97211  88.70925  87.91127
##            Aug       Sep       Oct       Nov       Dec
## 2022                               119.36102 118.06743
## 2023  93.61892 104.84220 115.93496 114.50585 114.71910
## 2024  92.04224 102.54024 115.82103                    
## 
## $upper
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 2022                                                                        
## 2023 128.6992 127.6015 125.1672 121.5254 116.4790 118.6528 118.4535 121.9950
## 2024 134.7366 131.8047 130.0801 125.6757 121.6150 123.3592 122.5642 126.6963
##           Sep      Oct      Nov      Dec
## 2022                   141.0694 143.8888
## 2023 133.2195 144.3128 146.7058 148.3774
## 2024 137.1948 150.4758
auto_fit$mse
## NULL
#[1] NULL
#RSME is 3.602064 on average forecast values were 3.602064 away from the actual
#MPE is 11% percenatge of error is around 11%


#accuracy summary

accuracy(ets_forecast)
##                       ME     RMSE      MAE        MPE     MAPE      MASE
## Training set -0.03996037 3.310275 2.498396 -0.1249896 2.482994 0.5687464
##                    ACF1
## Training set 0.07788466
accuracy(naive_forecast)
##                     ME     RMSE      MAE        MPE     MAPE    MASE      ACF1
## Training set 0.2536275 6.305638 4.625101 0.02809525 4.635846 1.05288 0.2207691
accuracy(rwf_forecast)
##                     ME     RMSE      MAE        MPE     MAPE    MASE      ACF1
## Training set 0.2536275 6.305638 4.625101 0.02809525 4.635846 1.05288 0.2207691
accuracy(snaive_forecast)
##                    ME     RMSE      MAE      MPE     MAPE MASE      ACF1
## Training set 2.731406 5.829941 4.392811 2.525765 4.342721    1 0.6254568
accuracy(SSE_Simple)
##                      ME    RMSE      MAE        MPE     MAPE     MASE      ACF1
## Training set 0.04457017 6.29552 4.615979 -0.1819999 4.631386 1.050803 0.2185404
accuracy(auto_fit)
##                        ME     RMSE      MAE        MPE     MAPE      MASE
## Training set -0.004444951 3.602064 2.583457 -0.1100037 2.571757 0.5881103
##                     ACF1
## Training set -0.06809078
#I would choose ets_forecast as it has the lowest MAPE. I chose MAPE because the units are in Percent, it is extremely useful 
#when observations are large numbers and also It can be used to compare same or different techniques as units is in %.


#Best MAPE Forecast

accuracy(ets_forecast)
##                       ME     RMSE      MAE        MPE     MAPE      MASE
## Training set -0.03996037 3.310275 2.498396 -0.1249896 2.482994 0.5687464
##                    ACF1
## Training set 0.07788466
#worst MAPE Forecast

accuracy(naive_forecast)
##                     ME     RMSE      MAE        MPE     MAPE    MASE      ACF1
## Training set 0.2536275 6.305638 4.625101 0.02809525 4.635846 1.05288 0.2207691
#conclusion 

#based on my analysis the value of time series will be cyclical as it has been so far in 1 and 2 both years